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Abstract. This study presents miseellaneous properties of pseudo-factorials, 
which are numbers whose recurrence relation is a twisted form of that of usual 
factorials. These numbers are associated with special elliptic functions, most 
notably, a Dixonian and a Weierstrafi function, which parametrize the Fermat 
cubic curve and are relative to a hexagonal lattice. A continued fraction expan- 
sion of the ordinary generating function of pseudo-factorials, first discovered 
empirically, is established here. This article also provides a characterization of 
the associated orthogonal polynomials, which appear to form a new family of 
"elliptic polynomials" , as well as various other properties of pseudo-factorials, 
including a hexagonal lattice sum expression and elementary congruences. 



1. Pseudo- factorials 

Start from the innocuous looking recurrence, 



which determines a sequence of integers with initial values 

(2) 1, -1, -2, 2, 16, -40, -320, 1040, 12160, -52480, -742400. 

These numbers will be called pseudo-factorials, since the omission of the sign alter- 
nation (—1)"+^ in Equation ([!]) determines the sequence of factorials, 0!, 1!, 2!, . . . 
At the suggestion of one of us (Bacher, October 2004), they have been included as 
Sequence A098777 in the On-line Encyclopedia of Integer Sequences |26j . which is 
brilliantly maintained by Sloane and a gang of dedicated volunteers. The purpose 
of this note is to show that these numbers, though not classical, have a host of 
interesting properties. 

The exponential generating function, 

(3) 1 + ^ + 2^ 

ri>0 

is fundamental to our treatment. (The fact that the absolute values |q;„| are domi- 
nated by factorials implies that f{z) is analytic at least in \z\ < 1.) We first 
elucidate the relation between f{z) and elliptic functions of a kind introduced by 
Alfred Cardew Dixon in 1890 {vide [7 ), then show the reduction to the more 
common WeierstraB framework of p- functions: this forms the subject of Sections [2] 
and [3] A simple consequence of the elliptic connections, worked out in Section ^ 
is an expression of pseudo- factorials as sums over a hexagonal lattice. Next, we 
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establish a continued fraction of a new type, relative to the ordinary generating 
Junction of pseudo-factorials, 

F{z) ^ a„z" = l-z-2z^ + 2z^ + 16z^ - 40z^ , 

n>0 

to be taken in the sense of formal power series since its radius of convergence is 0. 
To wit: 



(4) F{z) ^ 



1 



ra>0 



l + Z 



3 • P . z2 



1 - z- 



22. 



3 • 32 . z2 

1 + 3z + 



l-3z 



42.^2 



where the denominators are successively +1, — 1, +3, — 3, +5, — 5, . . ., and the nu- 
merators are 3 • I2, 22, 3 • 32, 42, 3 • 52, 62, . . .. Such a repetitive pattern of order 2 in 
what is known as a Jacobi fraction is somewhat unusual. The fraction Q was first 
discovered experimentally — proving it in Sections [5}]7] is the central theme of the 
present study. 

Finally, the convergents of the continued fraction Q can be made explicit, via 
generating functions: this is conducive to what seems to be a new class of "elliptic 
polynomials" in the sense of Lomont-Brillhart [T7j; see Section [sj In Section [9] we 
then draw several consequences of the previous developments, in the form of Hankel 
determinant evaluations and elementary congruence properties of pseudo-factorials. 

Why? As kindly suggested by a referee and by the editor of The Ramanujan 
Journal, we offer a few comments relating to the motivations behind studies, such 
as the present one. 

One primary motivation is to gain a better understanding of a still mysterious 
class of continued fractions, the ones that have coefficients that are polynomial (or 
rational) functions in the depth n. This highly general but somewhat impenetrable 
class of special functions seems to have been first noticed by PoUaczek in [22 . It 
is piquant to note that Apery's continued fraction |30j . which first led to an irra- 
tionality proof of Ricmann's C(3), has cubic denominators and sextic denominators, 
a feature also shared by the Conrad continued fractions [5J |B] relative to the Dixo- 
nian functions sm, cm, themselves closely related to pseudo-factorials (Section |2]). 
Our main continued fraction (Theorem |5]) is PoUaczek in a mildly extended sense: 
it has polynomial coefficients, but modulated by an unsual periodicity of order 2. 

On another register, the relation of elliptic functions to continued fractions is 
an old subject, going back at least to Eisenstein, Stieljes, Rogers, and Ramanujan. 
Whereas the q-series and -i^-function aspects are not immediately relevant to our 
discussion, we may observe that the continued fractions relative to the Jacobian 
functions sn, cn, dn have arithmetic content: together with the Hankel determinant 
product evaluations that they imply, they have been used by Milne in his penetrat- 
ing study [2D] of representations of integers as sums of 4n2 or 4n(n + 1) square or 
triangular numbers. In this spirit, we offer some new Hankel determinant evalua- 
tions in Section |9] but will leave to others the task of determining whether they are 
of some arithmetic interest. 
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The class of orthogonal polynomials associated to any continued fraction, whose 
coefRcients are rational-in-the-depth, is also of interest. Many instances of low 
degree have been categorized by Chihara in [3]. "Elliptic polynomials", that is, 
orthogonal polynomials associated with series expansions of elliptic functions have 
a tradition that goes back at least to Carlitz [2] and they form the topic of the 
entire book by Lomont and Brillhart fVP . The discovery of a new class (Section [s]) 
is of interest in this context. 

Finally, the pseudo-factorials are tightly coupled with ^T{^)^, as seen from our 
discussion of lattice sums in Section |4] In this context, it is well worth noting a 
spectacular arithmetic continued fraction for r(|)'^, recently obtained by Tanguy 
Rivoal [23j, which is of the sextic-duodecimic type(!). This and previous observa- 
tions reflect the fact that our understanding of an orbit of questions, surrounding 
Pollaczek continued fractions, elliptic functions, elliptic polynomials, and diophan- 
tine approximation properties, is still fragmentary; but they also suggest that "hid- 
den" structures are yet to be discovered in this area (see also our brief conclusion 
in Section 10). 



2. Elliptic connections: Dixonian functions 

This section serves to establish the first connection between pseudo-factorials 
and elliptic functions. The starting point is the exponential generating function 
defined by ([s]); it satisfies a functional equation, 

(5) f{z) = -f{-z)\ 

which directly translates the defining recurrence ([T|). 

To make f{z) explicit, take the functional equation ([5]) and differentiate once, 
so that 

(6) /"(z) = 2/(-2)/'(-z), implying f" {z) ^ -2^^^^) f{zf, 

since, by ([s]) again, one has f'{—z) = —f{z)^ and f{—z) = /'(z). In order to 
solve the nonlinear differential equation, "cleverly" multiply by ^J— f'{z) to get 

nz)^^nz) = 2f{zfnz), 

which is integrated to give 

(7) - I {-f{z)f' - Im' Ik, 

with K a yet unspecified constant. Equivalently, one has 

^ ' {K^f{zfY/^ 
which upon one more integration gives 

(9) / ... ...■.^./■. = 



(A'-u;3)2/3 

where use has been made of the initial condition /(O) = 1. The constant K is finally 
identified by means of a second order expansion (with /'(O) — —1, /"(O) = —2), to 
the effect that one must have K — 2. (The computations parallel those of [B].) 
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In view of our subsequent treatment, it is convenient to standardize A linear 
change of variables yields 



2-1/3 



(10) / = 



/2-V3/(.) (1-2/3)2/3 

where we have taken into account that K = 2. Throughout this study, a funda- 
mental constant is tts (a period of the function sm defined below in (14)), 

^''^ ^^'-='1 (w)^^^U'3j = r(iy = 2;^n3 

where the evaluation results from the classical Eulerian integral [5^ §12.4]: 
(12, B,..„:= -0'-<i<= r^^T^^"-™ 



(l + u)"+/3+i T{a + (3)' 

(Numericalljj^ we find 773 = 5.2999162508.) A simple computation shows that 

dy ^ l^fl l\ 



(l-y3)2/3 6 V3'3; 6' 
so that we can write 

The function f{z) can now be expressed in terms of specific elliptic functions 
introduced by A. C. Dixon in his original memoir f7^. Define the function sm(z) by 
the equation 



(14) 



(l_y3)2/3 



Then, a comparison of (13 1 and (14) permits us to identify f{z). 

Theorem 1. The exponential generating function of pseudo-factorials satisfies 

(15) /(z) = 2i/3sm(|-2i/3,), 



where the Dixonian elliptic function sia(z) is as in (14) and 1:3 is the constant (111. 

We can offer a few comments regarding Dixonian functions. There is actually 
a pair of "higher-order trigonometric" functions, sm and cm, where sm and cm 
are evocative of a sine and a cosine function, respectively. Their properties can be 
developed from first principles, as done by Dixon followed by Conrad-Flajolet [5j 
|5J|7], starting with the differential system, 

(16) sm'(z) = cm(z)2, cm'(z) = — sm(z)2, 

and initial conditions sm(0) = 0, cm(0) — 1. (See also the works of Lundberg [TS] 
and the recent developments by Lindquist and Peetre |15} 116] for a yet more general 
approach.) For the record, we note that 

4 7 10 3 6 9 

sm(z) = z-4|^ + 160^-20800^-- • • , cm(2) = l-2|^+4o|^-3680^ + - • • , 



^ The notation x = ^, with § a decimal fraction, indicates that ^ is an approximation of x to 
the last digit stated. 
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whose coefhcients (1, —4, 160, . . . and 1, —2,40, . . .) are respectively A104133 and 
A104134 of Sloane's Encyclopedia. The works of Conrad and Flajolet |H1 IE] pro- 
vide continued fraction expansions for the ordinary generating function of these 
coefficients, but these are then relative to the expansions of sm, cm at 0, and not 
at 713/6, as in (15 1. Finally, we observe that, by the calculation ([t]) and by (16 1, we 



have the identity 

sm(z)^ -I- cm(z)^ = 1, 

so that the pair (sm(z), cni(z)) parametrizes the Fermat cubic F3 de&ied by the 
equation + F'^ = 1, which is of genus 1. 



3. Elliptic connections: WeierstraC forms 

It is a priori possible to reduce the exponential generating function f{z) of 
pseudo-factorials to any of the several canonical forms of elliptic functions. Here, 
we show, by elementary calculations similar to the ones of the previous section, 
how to arrive directly at an expression involving the WeierstraB function p. We 
recall that this function p{z) = p{z; g2, gs) is classically defined by the nonlinear 
differential equation ^32J 

(17) p'{zf = 4p{zf - g2p{z)- gs, 

together with the initial condition p{z) as z ^ 0. By design, the pair 

(p(z), p'{z)) parametrizes the elliptic curve — 4X^ — g2X — g^, with invariants 

52,53- 

The starting point is the fundamental relation ([5|, namely, f'{z) = —f{—z)^. 
We first claim the identity 

(18) f{zf + f{-zf^2. 

The proof is obtained by verifying that the derivative of the left-hand side is 0, 
(/(z)3 + /(-z)3)' ^ 3/(z) V'(^) - 3/(-z)2/'(-^) = 0, 

(the final reduction uses ([s])), combined with the initial condition /(O) = 1. Next, 
we have 

(19) (-/(z)/(-z))' = /(-z)3-/(z)3, 

again by way of 
Now set 

giz) := -/(z)/(-z) = -1 + 3z2 - 3z^ + 3z^ - ^z« + ^z^" - ^z^^ + ■■■. 

The basic elliptic connection is provided by the differential relation 

(20) g'(z)2 = 4g(z)3 + 4. 



which is clearly of the WeierstraB type (17 1. To see it, it suffices to square the two 
sides of the identities (18 1 and (19 1, then compare the outcomes. Equation (20) 
then shows that g{z) is closely related to the elliptic curve £ defined by 

(21) AX^ + 4 . 
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The curve £ contains six integral pointfj^ (including the point at infinity correspond- 
ing to the identity element of the underlying group) forming a cyclic group of order 
six. The non-trivial points of this group are: (—1,0) (of order 2), (0,±2) (of order 
3) and (2, ±6) (of order 6). Since {g{0),g'{0)) = (—1,0), the series g{z) represents 
the expansion of the WeierstraB function p around the unique real 2-torsion point 
(-1,0) of£. 

The following result recovers f{z) from g(z) = —f(z)f{—z) and constitutes the 
main result of this section. 

Theorem 2. Let p{z) :— p(z; 0, —4) be the Weierstrafl function with invariants 
52 = and gs = —4 and smallest positive real perio^ 



(22) 



6,..,32-v>.?:^^rfl 

TT V 3 



with 11^ as in (111. The exponential generating function of pseudo-factorials satisfies 

-p'(z-|-3r) -2iV3 



(23) 



fiiViz 



2iV3p(z + 3r) 



Proof. We first establish the expression (22 1 of the real period of p(z), only making 
use of the most basic properties of elliptic functions (32 Ch. XX]. Let us denote 
temporarily the real half-period by zo. By general properties of elliptic functiontj^ 
we have p'(ra7) = 0, while p{w) is a real root of Aw^ -1-4 = 0; that is, p{ti7) = — 1. 
Thus, since p is the inverse of an elliptic integral, we must have 



dw 



6 



B 



1 

3' 



dw 



6 



B 



1 

1 1 

3 " 



1 



dw 



6 



-w'- 



and u 



(respectively) in the last 



as shown by the changes of variables t 
two integrals of the first line, followed by Eulerian Beta function evaluations (12 1. 
We henceforth denote the quantity w by 3r. 

The function f{i\/3z) is determined by /(O) = 1 and by the functional equation 
deduced from ([5|: 

(24) J_jL/(,V3,) = _/(„,V3z)^ 

iV3 dz 



We proceed to verify (23). Since p'(3r) = and p(3r) = —1, we first have 

-p'(3r) - 2^3 



2iV3p{3r) 



/(O) 



^ These are the only rational points of the curve £, since it is known that the Mordell curve 

= -(- 1 has six rational points; see, for instance, the SIMATH tables that are accessible at 
tnt . math .metro-u . ac . jp/siniath/MORDELL/MDRDELL+. 

■^The reader should be warned that the six numbers it6r, ±6re^^'" are not the shortest non- 
zero elements of the period lattice for p. They generate a sublattice of index 3 in the period 
lattice of p whose shortest elements are given by the six numbers it2j\/3r, i2iy/3re^^'^^^; see 
also Section |4] and Figure [l] 

^We have = since p' is odd; hence must be a root of the third-degree polynomial 

associated with p; here, p{zu) = 
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In view of ( 24 1 , the equality ( 23 1 then reduces to proving the identity 

2 



(25) 



1 -p'{z + Sr)-2iV3 
iVS \ 2^^/3p(z + 3r) 



-p'{-z + 3r) -2iVS 
2i%/3p(-z + 3r) 



= 0. 



Since p(z) is an even function, which is 6r— periodic, the left-hand side of (25) 
can be put under the rational form A/B, where the numerator A involves p, p', p" 
at z + 3r. Writing p, p', p" for p(z-f 3r), p'(2;-f 3r), p"(z-f 3r), it remains to verify 
that 

A = -2p"p + 2(p' + 2iV3)p' + p'^ - AlVSp' - 12 
vanishes identically. Derivation of the differential equation p'^ = ip^ + 4 for p 
yields p" = 6p^ and we indeed obtain 

A=-12p3-|-3p'2-12 = 0, 

which concludes the proof of the statement. □ 



4. Lattice-sum expressions for the pseudo-factorials 

From the Dixonian as well as the WeierstraB connections discussed in the previ- 
ous section, expressions of the q;„ as lattice sums can be developed. 

Theorem 3. The pseudo-factorials are expressible as lattice sums involving a 
twelfth root of unity: with p — 2r\/3 and r as in Eq. (22 I, one has, for any n > 2: 



(26) ""--3^TT 2^ 77\ IX ' ^ 



[(A-^)C+(A^-i)C- 



e 



27r/6 



The formula ( |26[ ) implies explicit asymptotics for (a„), namely. 



z — (p + tj) p + cj (p + w)2 



{2i^y. ' ' {2iy+iy. 

once only the relevant dominant poles are retained. This explains in particular the 
regular sign pattern "H h" observed in (|2|. 

Proof. Let G{z) be an elliptic (i.e., meromorphic, doubly periodic) function that 
has only simple poles. Let A = Zuji Zuj2 be its lattice of periods and let So be 
the set of poles contained in a fundamental domain of A. Then, if ^ Sq, one has 

(27) G{z)^G{0) + zG'{0) + J2T.^P ^ ^ 

where represents the residue of G{z) at z — s. Theorem [2] and the formula (27) 
show that it suffices, up to an affine transformation, to work out the singular 
structure of 

H{z) = -^'(^) - 2iV3 
2i^/3p{z) 

in order to deduce the partial fraction expansion of /(z), from which the lattice 
sum (26) expression will result. 

The function p(z) = p(z; 0, —4) has lattice of periods 

(28) {Ze''^/^ © Ze"'"'/^)2\/3r 
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Figure 1. Left: the "primary" lattice of periods [tliick lines] of f{z), 
where a fundamental domain is obtained by the union of two adjacent 
equilateral triangles (one grey, one white); the poles (three per funda- 
mental domain) are represented by small discs and form a "secondary" 
hexagonal lattice [thin lines]. Right: a diagram showing the three poles 
(circled) of a fundamental domain of f{z) around the origin. 



and H{z) has simple poles at z = and z = ±2r. Given the series expansion 
p{z) = z~'^ + 0{z'^) of p around 0, wc find the residue of H{z) at to be 

(29) Res(i7(z);z = 0) = ^. 

By the discussion following ( [2T| ), we also have (p(±2r), p'(2r)) = (0, ±2). Moreover, 
we have (p(3r), p'(3r)) = (-1,0) and (p'(3r), p"(3r)) = (p'(3r), 6p(3r)2) = (0,6), 
hence p'(2r) = —2 and p'(4r) p'(— 2r) ~ 2. The expansion of p{z) around 2r is 
thus 

p(z + 2r) = -2z + 2z* - ^z^ + ^z^" + ■■■ , 

from which a simple computation provides the residue of H(z) at 2r, and, similarly, 
the residue at — 2r: 

(30) Res (H(z) ; z = 2r) = ^ + ^ Res (71 (z) ; z = -2r) = ~^ + ^ 

6 6 



Now, by (28), (29), and (30 1, the singular structure of H(z) is entirely known. 
Since /(izyo) = H{z + 3r), an afhne transformation (composed of a translation, 
a rotation, and a dilation) provides the singular structure of /(z) itself — we ab- 
breviate the discussion, which is routine. As represented in Figure [T] the lattice of 
periods A of /(z) is a hexagonal lattice with generators 3r(— l±i-\/3); we may call it 
the "primary" hexagonal lattice. There arc three poles of /(z) in the fundamental 
domain, at — 3r and zti\/3r, which, upon translation by A, generate a "secondary" 
hexagonal lattice. The residue of /(z) (deduced from (|29| and ([30|), is then found 
to be of the form (8^+*/^ at a point of the form 2rV3[{X - l/2)( + {fi - 1/2)C"^], 
which corresponds to a 3-colouring of the secondary hexagonal lattice (since is a 
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third root of unity) . The corresponding partial fraction decomposition ( 27 1 results 
for f{z). Finally, the faci[^that, for n>2, 

r ni 1 1 z \ 1 

] + - + 



— a a a? I 

yields the stated lattice- sum formula for a„. □ 



The sum ( 26 1 establishes the pseudo-factorials as a two-dimensional analogue 
of Bernoulli and Euler numbers, one that is relative to the hexagonal lattice. It 
might be of interest to investigate systematically continued fractions relative to 
other hexagonal lattice sums. It is worthy of note that arithmetic properties of 
analogous "lemniscatic" sums, relative to the square lattice, have been studied by 
Hurwitz [Hin]. 

5. The Stieltjes— Rogers addition theorem and continued fractions 

This section serves to introduce the basic technology needed to develop an ex- 
plicit continued fraction expansion from an addition theorem of a suitable form. An 
experimental approach specialized to pseudo-factorials follows, in Section|6] We can 
then reap the crop in Section [7] and finally prove our main continued fraction result 
(Theorem [5|. 

Stieltjes and Rogers independently discovered that the continued fraction ex- 
pansion of an ordinary generating function, <I>(z) = ^(/)„z", is closely related to 
addition formulae satisfied by the corresponding exponential generating function, 
= X] '^n-^"/"-'- First, a definition. 

Definition 1. Let 4>{z) = 1-1- X]n>i 4'nz"' /nl be a formal power series. It is said 
to satisfy an addition formula of the Stieltjes-Rogers type if there exist nonzero 
constants uji, UJ2, ■ ■ ■ and formal power series ipq(z), ipi(z), . . ., such that 

(31) (l){x + y) = 'fio{x)^o{y) + '^iV'i{x)(pi{y) + uj2V2{x)(p2{y) H — , 
where (pr{z) has valuation r and is normalized by (pr{z) = {z^ /r\) + 0{z'^'^^). 



In (31 1, the valuation condition on Lpr is essential, the normalization (firiz) ~ 
z'^ /r\ being a mere convenience for what follows. 

The addition formula gives rise to an algorithm for computing the (pi{z), know- 
ing (t>{z), either symbolically or via some series expansion. First, setting y — m 
the addition formula shows that (pi^z) must be equal to (po{z) (this makes use of the 
normalization 1^90(0) = !)• Next, assume that the functions (/jq, ■ • ■ , </'^-i and the 
coefficients lui, . . . jW^-i have already been determined. Then, by differentiating £ 



times the addition formula pl| with respect to y, then setting y = 0, one finds 

9i<P{x) - ^ujjipj{x) [dl(Pj{y)]y^g = uje(pe{x). 
j=o 

Upon comparing the coefficient of x^ in the Taylor expansions of both sides, we 
see that at most one nonzero coefficient ajg can satisfy the equation, given the 
normalization (pi{x) — /£l + 0{x^^^). Once this choice has been fixed, then (pe{x) 



^ Notation: we let [2"]/ represent the coefficient of z" in the formal power series or analytic 
function /. 
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is uniquely determined as a linear combination of the previous ipj{z), together with 
derivatives of 0, and the process can continue. This construction also determines 
the broad class of functions in which the >Pi{z) live: each <pe{z) belongs to the vector 
space generated over C by the first derivatives 0, (j)'^^^ . . . , (/'^^■' of the function (j){z). 
This algorithm, albeit suboptimal from a computational point of view, permits us 
to experiment with addition formulae relative to the generating function of any 
given number sequence (0n), a fact that will prove especially valuable in Section |6] 
Addition formulae of the Stieltjes-Rogers type are logically equivalent to contin- 
ued fraction expansions as expressed by the following central theorem originally due 
to Stieltjes [29] and Rogers [21; see also Perron [21J p. 133] and Wall [31l p. 204]. 

Theorem 4 (Stieltjes-Rogers). Let the exponential generating function (f>{z) = 



1 + X]ri>i 4'nz"'/n\ satisfy an addition formula of the form (31 1. Then, the corre- 
sponding ordinary generating function ^{z) = 1 + X]n>i 4'nz"' admits a Jacobi-type 
continued fraction expansioT^ 

(32) <^>iz) = ^ ^ , 

aiz 

1 - cqz ^ 

1 — CiZ 



where the coefficients are determined by 
a, = ^ (J>1), c 

There, (pj^k — k\[z'']^pj{z) , (p~i^k = 0, and cjq = 1 

Stieltjes am 
(z) = sec(2:) 



= — ^ (j > 1), Cj = (y5j,j + l - (j>0). 

^j-i 



As an illustration, following Stieltjes and Rogers, we examine the case of 

1 



cos(z) ' 

Each (pk{z) (provided it exists) must then a priori be of the form sec(z)Pfe(tan(z)), 
where Pk is a polynomial satisfying deg Pk = k. In a simple case like this, classical 
trigonometric identities yield 

sec(x + y)= ^ \ = Y] sec(a;) tan(a;)'' • sec{y) tan(y)'=, 

cos(a;) cos{y) — sm(a;) sm[y) 

which, in normalized form, becomes 

(33) secix + y) = Y^ikir (secix)'-^) ■ (sec(y)^)' 

fc>0 ^ • / \ 



With En = n![z"] sec{z) an Euler number. Theorem |4] then provides the continued 
fraction expansion: 



1 



z 



n>0 I 



22. z2 

1 



^Such an expansion is diversely known as a Jacobi fraction, a J-fraction, or an associated 
continued fraction. 
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where the coefficients 1^ ,2^ ,3^ , . . ., are obtained here as quotients of consecutive 
squared factorials. The absence of linear terms reflects the parity of sec(z). 

6. Experimental determination of the addition formula for f{z) 

Section |5] has shown that, in order to approach the construction of a continued 
fraction expansion for pseudo- factorials, we need to develop a suitably constrained 
addition formula for their exponential generating function f{z), which is elliptic. 
We proceed here in an experimental manner in order to infer the likely shape of such 
an addition formulsj^ Once this had been done, the proof of our main continued 
fraction reduces to purely mechanical verifications to be carried out in the next 
section. 

A first idea that comes to mind is to look for an addition formula of a kind 



similar to the secant case ( 33 1 , namely, 

(34) <j}{x + y)^ <j)ix)(biy)^{c7ix) ■ a{y)), 

for some function (or power series) ^'(w). However, all such formulae can only 
arise from a class of special functions that comprises five parametrized subclasses, 
of which prototypes are 

sec(z), e-'-\ e^'/2 



1-z' ' ' 

(This is a rephrasing of a classification of orthogonal polynomial systems due to 
Meixner [3j 19 .) Obviously, elliptic functions are not amongst this group. 

Another source of inspiration is a continued fraction relative to elliptic functions, 
which is also due to Stieltjes and Rogers. With sn,cn,dn the Jacobian elliptic 
functions, as classically defined, where we leave the modulus k implicit, we have 

cn a; cn y — sn a; sn J/ dn a; dn ?/ 

(35 cn(a; + y) = — ^ ^ , 

1 — fc^ sn^ X sn^ y 

see for instance [321 p. 497]. This can be put into an equivalent Stieltjes-Rogers 

form (up to normalization) , namely, 

(36) 

cn{x+y) = cnxcn?/ (l — sn xdn a; sn ydn y + fc^ sn a: cn^ x dna:sn?/cn2 ydn y + • • • ) 
corresponding to the continued fraction expansion 

^ cn„ z" = l2T^2 ' 

n>0 I 



where cn„ n![z"]cn(z). 

We now turn to the continued fraction expansion relative to pseudo-factorials 
which, by Theorem |4] involves determining the right addition formula for f{z). 
Based on experiments under the Maple system as well as on induction from the 



^This section is not, strictly speaking, necessary. It could have been replaced by the shorter 
but somewhat obscure formulation: "Crystal ball gazing revealed to us the addition formula BOl." 



12 



ROLAND BACKER AND PHILIPPE FLAJOLET 



secant (33 1 and Jacobian (35 1 cases, we started searching for an addition formula 



of the form 

(37) fix + y) = fix)f{y)^{a{x)aiy)) + h{x)h{y)E{Tix)T{y)), 

for some power series S, /i, a, r, with (at least) ^'o = So = 1, and (t{x), t{x), h(x) 
all being 0{x). Since some binary pattern is present in the continued fraction, 
it is natural further to suppose that (j{x) = 0{x^), t{x) — 0{x^), which then 
corresponds to an "odd-even" addition formula, 
(38) 

f{x + y) = f{x)f{y) + h{x)h{y) + 'i'if{x)f{y)a{x)a{y)+E,h{x)h{y)T{x)T{y) + - ■ ■ , 



where — [w"^]"^{w) and = [w"^]'E.{w). In the notations of (|31J), we thus 
hope for an addition formula of the form 

(39) (p2j{x) cc f{x)a{xy , ip2j+i{x) cc h{x)T{xy~'^ , 

where a{x) oc b{x) means that the ratio a(x)/b{x) is a constant. 



The pleasant feature of the conjectured expansions ( 37 1 and ( 39 1 is that their 
plausibility can be effectively tested. Indeed, from the previous section, we have 
available an algorithm that can determine the (unique) y^j{x) corresponding to 
(j){z) = f{z), this to any desired precision. It then suffices to check that 

ip2{x) ip4{x) (p3{x) (p5{x) 

— --^ (X — --^ cx • • • , and — — — cx — --— cx • • • . 

(po{x) (p2{x} (pi{x} ip3{x) 

Verification of these relations for about a dozen of the ipj and till orders in the 
range 50-100 convinces us that we are on the right tracks. 

In fact, we found experimentally that a{z) = t(z) up to 0(z^°), the series 
starting as 

a{z) - 3 fz^ - + z6 - -z8 + -zi" - . 
V 7 7 



Also, the function h{z) that appears in (38 1 must be proportional to </?i(z) of the 



addition formula, whose expansion starts as (pi{z) = z — z^ + ^z"* — • • • . That 
function </?i(z) must itself, on general grounds, be a linear combination of /(z) and 
/'(z) without constant term, so that 

</Ji(z) = (/(z) + /'(z)) and h{z) (x (pi{z). 

Finally, assuming the observed law of the coefficients in the continued fraction to 
hold forever, we can deduce the only possible shape of the S and ^' functions. The 
function a is then inferred on the basis of the fact that (p2{z) cx ipo(z)a{z) must be a 
linear combination of /, /', /". All in all, every ingredient of an addition formula of 



type (37 1 is in place, and we are eventually led to conjecturing an addition formula 



for /(z) 



f{x + y) 



f{x)f{y) - lh{x)h{y) 



(40) <( l^la{xMy) 

h{z) = /(z) + /'(z), ct{z) - l-/(z)/(-z), 

which we shall establish in the next section. 
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7. Proof of the continued fraction expansion 

At this stage, we know that estabhshing the continued fraction Q relative to 
the ordinary generating function F{z) of pseudo-factorials reduces to deciding the 



validity of the conjectured addition formula (40 1 for the exponential generating 
function /(z). The proof we propose is a computer-assisted verification. As we 
shall explain, it only involves routine algebraic manipulationsj^ namely, rational 
function operations, normalizations, substitutions, as well as multivariate polyno- 
mial divisions. The calculations were performed using the Maple computer algebra 
engine (version 11). Without any attempt at optimization (we purposely wanted 
our program to rely solely on the most basic algebraic operations), the mechanical 
verification reduces to the mere execution of a few billion machine instructions — 
currently, just a few seconds of elapsed time. 



Proposition 1. The function f(z) satisfies the following addition formula: 
(41) f{x + y) 



f{x)f{y)~\{f{x) + nx))u{v) + rm 



l-\{l-f{x)f{-x)){l-f{y)f{-y)) 

Proof. We can a priori appeal to either the Dixon or the Weierstrafi framework, 
and we have opted for the latter. The WeierstraB p-function, p{z) = p(z; 0, —4) 
satisfies the two algebraic relations 

p'(z)2 = 4p(z)3+4 (DEF) 

The first one (DEF) is the basic differential equation, which serves as definition 
of p; the second one (ADD) is the familiar addition theorem of elliptic function 
theory [311 P- 441]. Both are "known" to Maple; both can be viewed as determin- 
istic rewrite rules permitting one to expand and simplify expressions involving p. 

Let 6r be the fundamental constant (real period) of Section [s] We know that 
(p(3r), p'(3r)) = (—1,0). The addition rule (ADD) combined with the expression 
of f{z) stated in Theorem [2] Equation ( p3| , then mechanically expresses f{z) as 
a rational fraction in p{Z) and p'{Z), where Z = z/iV^. Similar expressions are 
obtained for /'(z) (by standard derivation rules combined with partial reductions 
by (DEF)) and f{~z) (since p is an even function, while p' is odd). In this way, 
one automatically obtains rational forms in p, p' for 

fix), fix), fiy), f'{y), fi~x), fi~y), f{x + y), 

where the last one necessitates a substitution z i— s- x + y, followed by an application 
of the (ADD) rule. 

Let now T> be the difference between the left-hand side and the right-hand side 



of the relation to be proved. Equation (411. By the process described above, V 



becomes a rational function, with coefficients in Q{i^/3), in the four quantities 
Xi,Yi,X2,Y2, where Xi — p{x/i^/3), Yi = p'{x/i^/3), and similarly for X2,Y2, 
with y replacing x. The (rather large) rational fraction normalizes to the form 



^The validity of intermediate steps is, in addition, easily cross-checked by means of Taylor 
series expansions. 
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T) = A/B, with numerator A and denominator B involving, respectively, 2388 and 
1256 monomials. One can then operate with the rule (DEF), instantiated as 

^ iXf + 4, Yi 4X| + 4, 

the corresponding reductions being simply effected by multivariate polynomial di- 
visions. When this is done, we find that A reduces to 0, while B is reduced to a 
nonzero polynomial, which is of degree 1 in Yi , Y2 and of degree 8 in Xi , X2 ■ The 
verification of the addition formula for /(z) is thereby completed. □ 

A direct application of Theorem [4] to the addition formula expressed by Propo- 
sition [T] gives rise to our main continued fraction. 

Theorem 5. The ordinary generating function of the pseudo-factorials satisfies 
(42) F{z) ^ YI = JTl^^ ' 

n>0 l + Z-\ 



2' ■ z 

1 -z+ 



3 • 32 . z2 

l + 3z+ 



42.^2 
l-3z+ 



where the coefficients are, with the notations of the Jacobi form (32); 
(43) c, = i-iy-' (^j + , a, ^ -f{2 - (-1)^) 



Proof. We make use of the addition formula (41 ), which, taken under the form (40 1, 
yields 

00 00 

n=0 n=0 

with the LiJn and ipn determined by 

^2„¥'2n(x)<^2„(2;) - +i-"f{x)f{y)a{xra{yr 

i^27i+W2n+i{x)'P2n+i{v) = -3"""^ft,(a;)ft,(y)a-(x)"(T(j/)"' . 

We have, with the notations of (40 1, f{z) = 1 — z + 0{z'^), as well as cr{z) = 



1 - f{z)f{-z) = 3z2 -I- 0(z4) and h{z) = f{z) + f{z) = -3z -I- 0{z^). The required 
normalization of a Stieltjes-Rogers addition formula, iy9„(z) = z"^ /n\ -\- 0(z"+^), 
combined with the low-order expansions of f (z) , h{z) , a (z) , gives us 

^2„ = 3"(2n)!2, c^2„+i = -3"+i(2n + 1)!^, 

as well as 

_2n 2n+l „2n+l 

We thus have 

^^ = -3, '^ = -22, '^ = -3.32, '^ = -42,..., 

LUq UJl UO2 

and, for Cj — ^j-i.j- 

Co = -1, ci = 1, C2 = -3, C3 = 3, C4 = -5, C5 = 5, • • • . 
By Theorem |4] the last two formulae conclude the proof of ( 42 1 . □ 
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8. A family of orthogonal polynomials 

Our goal in this section consists in finding an explicit form for the polynomials 
that appear in the convergents of the main continued fraction of Theorem [s] this 
by way of their exponential generating function. We focus our attention on the 
denominator polynomials, precisely, on their reciprocals, which form a family of 
formally orthogonal polynomials that appears to be new. 

We start by specializing to the continued fraction under consideration ( 42 1 some 
well-known algebraic properties found in [211 131] that hold for an arbitrary Jacobi 
fraction ( 32 1 . The convergents of ( 42 1 are obtained by truncating the infinite frac- 
tion before a numerator. In this way, a collection of rational fractions Pk{z)/Qk{z) 
of increasing degrees is obtained, 

1 1 - z l + 2z + z^ 1 - z 22 z2 - 30 



1 ' l + z' l-f 2 z2 ' 1 3 z -f 6 z2 10 z3 ' IH- 24 z2 - 8 z3 24 z4 ' 

so that Qq = 1, Qi = 1 + Z; and so on. The denominator polynomials Qk satisfy a 
"three-term recurrence" relation, 

(44) Qk = (1 - Cfe_iz)Qfe_i - ak-iz'^Qk^2- 

(The Pk satisfy the same recurrence, but with initial conditions Pq = 0, Pi = 1.) 
The reciprocal polynomials defined by 

(45) qk{z) = z'^Qk Q 
then satisfy the recurrence 

(46) = (z - Ck^i)qk-i - afc-i'7fc-2, 9-i = 0, go = 1, 



with the a.j,Cj as in (43 1. On general grounds, they are formally orthogonal 
with respect to a (formal) measure whose moments coincide with the pseudo- 
factorials, (a„). In other words, they are orthogonal with respect to the bilinear 
form 

(/,5> = (/5>, with (z")=a„. 

Observe finally that, once the Qk are known, the Pk can somehow be regarded as 
known. Indeed, relative to (32 1, one has, in the sense of formal power series, 

Qk{z)F{z)-Pk{z) = 0{z^''), 

so that the coefficients of the Pk are expressible as a convolution of the two sequences 
[z"]Qkiz) and a„ = [z"]P(z). 

We have the following characterization. 

Theorem 6. Let T(z, t) be the exponential generating function of the reciprocal 
polynomials (qk) of ( [is] ), with coefficients (43 1; 

fc=0 

Consider the algebraic curve 

(47) 2 + 3< + 3i(l + <)?7-2(l-3t2 + 3t^)?7^ = 0, 
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which is of genus and is parametrized by 

1 (w^ + 3) w ^{w + 1) + 1) 

3 + 1 



(48) ^ = o , i ^ V = i 



and let r](t) be the branch that satisfies 77(0) — 1; 



^2 t^ t^ t^ f 

7j{t) = 1 + i + 2- + 10- + 24- + 280- + 400- + 12880- 



Define 



/ . so 2t(l + t) t^ t^ r 

(49) xit) ^vit? - Y^WTW - ^ + 2! - ^3! + 4! - 1««5! ' ^^^6! 
and introduce the fundamental elliptic integral 

du t'^ t^ 

(50 J(t):= = = t + 3- + 45- + 1215- +8505- - ■ 

^ ' io VI - 3u2 + 3^4 3! 5! 71 91 

Then, the generating function T satisfies 

(51) T(z, <) = 77(i) cosh (z J(i)) + xit) sinh (z J(t)) . 



Equation (51 1 was first arrived at by a combination of induction and of partly 
heuristic calculations, based on "guessing" intermediate differential equations as 
well as on Maple's symbolic integration capabilities. Rather than offering a heavy 



proof by successive transformations of the defining recurrence ( 46 1 , we have opted 



to present a computer-assisted verification of (51 1. In this way, we feel we save 
symbols, hence pages, hence trees. The price to be paid was only a few hours of 
interaction with the symbolic engine and (eventually!) a few seconds of computer 
processing time. As in the previous section, only well-specified totally-algorithmic 
steps are eventually used. Once more, there is no difficulty in checking intermediate 
steps against series expansions up to order 100 and beyond. 



Our proof of the identity (51 ) eventually boils down to exhibiting a fourth-order 



differential operator in t, with coefficients in C(z), that is satisfied by the difference 



between the two sides of (51 1. It is then sufficient to check that both sides satisfy 
the same initial conditions given by the coefficients of i° up to t^. 

The entire process relies on the holonomic framework pioneered by Zeilberger |33j . 
with supporting theorems to be found in works of Stanley, such as [2T and PF', 
Ch. 6]. Let K be a ground field, which we take here to be C(z), the field of rational 
fractions in z. (Throughout, we treat the quantity z as a parameter.) A formal 
power series of Kp], simply called "function", is holonomic (alternative names 
are differentiably finite, D-finite, (9-finite) if it satisfies a linear differential equa- 
tion with coefficients in the rational field IK(t). Equivalently, h is holonomic if the 
vector space over K(t) spanned by all the derivatives {dj h} is finite-dimensional. 
Holonomic functions are known to be closed under sum, product, differentiation, 
integration, and algebraic substitutions (i.e., substitutions of algebraic functions in 
place of variables). Finally, if h is holonomic, its sequence of coefficients ([z"]/i) 
satisfies a linear recurrence relation with coefficients in K(n). 

Clearly, a holonomic function is determined by a finite amount of information; 
namely, a defining differential equation supplemented by sufficiently many initial 
conditions. Given two holonomic functions A, B, one can then verify their conjec- 
tured identity as follows. 
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(i) Compute a differential equation, of order lu, say, that is satisfied by the 
difference A — B. 

(ii) Check the coincidence of the expansions of A and B up to terms of order 

o(r). 

By the finiteness of the underlying vector spaces, the process constitutes a valid 
proof of A — B, in "non-singular" cases at leasij^ In our context, it could be carried 
out comparatively easily, thanks to the powerful Gf un library developed by Salvy 
and Zimmermann [25 . 

Proof (Theorem^. In what follows, we use d = dt to represent the differential 
operator we denote by S = S„ the shift operator on infinite sequences such 
that S(u„) = Un^i. We let 11^ generically represent a polynomial of degree r, either 
in t (for differential operators) or in n (for difference operators), with coefficients 
in K. As indicated before, the quantity z is treated as a parameter. Our purpose is 



to prove A ^ B, where A is the left-hand side of (51 1 and B is the right-hand side: 



(52) A^B, with A := T(z,t), B ?7(t)cosh(zJ(<)) + x(^)sinh(zJ(i)) , 

and r]{t), x(^)j J{i) ^ defined in the statement. See Figure [2] for a summary of the 
main steps of our proof. 



The left-hand side (A). The parity inherent in the coefficients (43) suggests to 
introduce the subsequences r„ = g2n and s„ — q2n+i- The basic recurrence (46) 
then relates r„ to r„_i, s„_i and s„ to r„, s„_i; hence, by substitution, the fact that 
the vector (r„, s„) depends linearly on (j-n-i, s„_i) via a matrix, whose coefficients 
are polynomial in n (and the parameter z). By instantiating this last relation 
at n + 1 and n + 2, and using back substitution, there results that r„ and s„ satisfy 
explicit linear recurrences of order 2 with coefficients that are polynomial in n. The 
difference operators annihilating (r„) and (s„) are found in this way to be of the 
form 

Hi +n3 +n5 S° . 
Equivalently, the exponential generating functions 

are found to satisfy £[i?] = and 9Jl[S'] = 0, where £ and 971 are each of the form 

In fact, second-order operators £° and 3Jt° that appear to cancel R and S, respec- 
tively, can be guessed (roughly, by the method of indeterminate coefficients cleverly 
implemented in Maple's Gfun). The guesses can then be turned into full-fledged 
proofs by checking (with Maple's Ore_algebra, see [1]) the operator divisibility 
relations: £° | £ and 971° | 971. Once this is done, a differential operator that an- 
nihilates T can be obtained by making use of properties of holonomic functions 



operator is said to be "non-singular" if the lead polynomial of the associated recurrence 
(relating the coefficients hn := [t"\h{t), of a solution h{t) and having coefficients in K[n]) has 
no root in Z>o. In the case of a non-singular operator of order u), the number of needed initial 
conditions equals lo. (In the "singular" case, a higher, but still effectively computable, number 
may be needed.) 
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n4d-'Y + U3dY - z^Y = 0; see (54 1 



A = T{z,t) : Ued'^Y + U5dY + UiY = 0; see 21 in ([531 

Vit),x{t) ■■ Ued^Y + U5dY + U4Y = 0; r?, x algebraic (|47f, (|49l 

1 exp(±zJ(t)) I 
1 cosh, sinh (2 J(t)) J 

f B+ = r7(t)cosh(zjm) , 

<^ y. Hiaa^y + '-' + ngy = 0; closure algorithm (x) 

\ B- =x{t)smh{zJ{t)) ' . 5 ^ ; 

A — ^ - {B+ + B') ni2d'^Y + h UgY = 0; closure algorithm (+) 

Figure 2. The shape of the differential equations satisfied by quanti- 



ties intervening in the proof of the main equation (51 1 



(effective closure under sum) and a further round of simplification based on guess- 
ing. It is found in this way that the second-order operatoip"] 
(53) 

a = 4 (1 - 3*2 st^){{2t{t + i)c - 1) a" 

+ 4 {24:t^( + 30t^( - (18 + 6 C) - 12 1\ - (4 C - 9) t - 2 C) 

+ (48t*C + 72i\- (48-H8C^ - 18C) t'^ - + SC^ + 12) t + 9 + + 12(^) d° , 

with C 2 — 1/2, annihilates T{z,t). 

The right-hand side (B) . We can build up differential equations starting with the 
explicit expression of i? in ( 52 ) : see again Figure [2] for a summary. Given a quantity 
X that depends on z, we set = \X{z) + ^X{—z) and X~ = \X{z) — ^X{—z), 
defining its "odd" and "even" parts (in z), respectively. With B = T{z,t), we 
then consider B^ = T"*" = ?7(i) cosh(z J(i)) and B^ = = x(t) sinh(z J(t)), and 
proceed to construct the corresponding annihilators, and 05^. 

First, we observe that if P(t) is an arbitrary polynomial, then 
(54) Y(t) exp I z / | satisfies Pd^Y + -P'dY - z^Y ^ 0. 

This equation is invariant by z ^ — 2, so that it is also satisfied when the exponential 



in (54 1 is replaced by sinh, cosh. 

Next, the function ri{t), given by a cubic algebraic equation, is found to satisfy 
a second-order differential equation with coefficients that are of degree at most 6. 
The application of closure rules for products of holonomic functions then provides 
for = T+ a differential operator 03"'" that is of order 4, with coefficients of 
degree at most 12. 

We can then proceed to construct the annihilator of the odd part B^ = T^. 
It turns out that the algebraic function x(t) defined in (49) satisfies the same 



differential equation as rj{f), but with different initial conditions (x(0) = 1, x'(0) = 



0). There now results from this fact and the comments accompanying (54 1 that we 
can take 

<8 = 05+ = OS", 



^'^ Much to our surprise, Maple's symbolic integrator proposed a solution to 2I[/] = 0, which 
involved terms of the rough form exp {±zj{t)) and eventually led us to infer (|51|. 
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as annihilator of the right-hand side (B) of (51). 

The comparison. Finally, it remains to verify that A ^ B. The operator A is 
defined to annihilate A -|- B, where A and B are the vector spaces of solutions of 
2l[/] = and = 0, respectively. By construction, the difference A — B is such 
that A[A — B] = 0. The operator, obtained by holonomic closure under sums, is of 
type 

A = Ui2d^ + Uiid^ + Uiod^ + Ugd^ + Uad°. 

The associated recurrence operator is found to be of the form 114 + • • • + 114 S", 
with leading coefficient 

(n -t- 9)(n -I- 10)(n -I- ll)(n -|- 12)(9 -I- Az^), 

so that A is non- singular at 0. It thus suffices to verify that the expansions of 
A = T{z,t) and of the right-hand size B in ([sT]) coincide till terms of order 0{t'^), 



A,B =l + iz + l)t+ (z^ + 3z^ + 6z + l0)*- + {z^ + 2Az^ - 8z + 24) ^ + 0{t^), 



so as to complete the proof that A = B. Equation (51 1 is now established. □ 



Orthogonal polynomials attached to continued fraction expansions relative to 
elliptic functions, have been first studied by Carlitz and Al-Salam (see [T31 HI] for 
some more recent developments), and, as already mentioned, they form the subject 
of the monograph Elliptic Polynomials by Lomont and Brillhart |17j . The family 
made explicit by Theorem [6] does not appear to be captured by their classifica- 
tion and hence seems to be new. Remarkably, in connection with birth-and-death 
processes having cubic weights, Gilewicz et al. [5] have recently discovered another 
new family of orthogonal polynomials, related to the expansion of the Dixonian 
function sm taken at (as in [B^), rather than at the point tts/G that is needed here 
(cf Theorem [l]). 



9. Consequences of the continued fraction expansion 

The continued fraction of Theorem |5] has several interesting by-products that we 
now examine. These include an explicit evaluation of Hankel determinants, as well 
as elementary congruence properties of pseudo-factorials. 

Hankel determinants. It is well known that, generally, coefficients of a Jacobi 
fraction can be expressed as determinants. This fact is classically derived from 
Stieltjes's matrix version of the addition theorem [3TJ pp. 202-206]; it is equivalent 
to the LDU decomposition of the Gram matrix H, with entries hij — {z'^,z^) = 
(for i, j > 0), which is also known as the Hankel matrix of the sequence (z"); 
see for instance |12l §2.1]. Conversely, any known continued fraction yields an 
explicit Hankel determinant evaluation. Given this, an immediate consequence of 
Theorem [5] is the following. 

Corollary 1. Let m be a positive integer. The Hankel determinant of pseudo- 
factorials 

ao ai ■ ■ ■ am~i 
ai a2 ■ ■ ■ a-m 

Ctm-l CXm ■ ■ ■ 0.2m-2 
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Figure 3. The congruences (a„ mod M) for M = 2, . . . , 20 and n = 0, 
admits the closed form 



,25. 



H^r^' = n 



(m even) 



where the aj = — j (2 — (—1)-') are the continued fraction numerators of (43l 



Congruences. A cursory examination of the a„ suggests clear divisibility patterns; 
for instance, from (|2]), we immediately expect the a„ to be divisible by 10, for 
n large enough. Figure [3] tabulates arithmetic congruence properties of the a„ 
for small values of the modulus M and of the index n. The table obviously has 
much structure: the sequence (q!„) appears to be eventually modulo the numbers 
2, 4, 5, 8, 10, 11, 16, 17, 20; there are obvious periodically reproducing patterns, such 
as 1,2 (mod 3), 4,2 (mod 6), or the more recondite, and curiously repetitive, 

6,5,2,2,2,2,4,1,6,6,6,6,5,3,4,4,4,4,1,2,5,5,5,5,3,6,1,1,1,1,2,4,3,3,3,3 

of length 36 corresponding to modulus 7. We state here a simple consequence of 
our main continued fraction ( 42 1 in Theorem [5] 

Corollary 2. (i) The sequence (an) of pseudo-factorials is eventually periodic mod- 
ulo any integer M > 2. (ii) For each m > 2, the sequence (an) satisfies mod- 
ulo M ~ 3^™/^^™!^ a linear recurrence with constant coefficients that is of order at 
most m. 



Proof (sketch). The statement is an instance of the general fact that J-fraction 
expansions with integer coefficients automatically imply congruence properties [S]. 
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(i) The main continued fraction (42) representing the ordinary generating func- 
tion F{z) of pseudo- factorials has a factor of at its numerator of rank {M + 1). 
In particular, the contributions induced by the stages {M + 1), {M + 2), and so on, 
of this continued fraction are zero modulo Af^. In other words, F{z) is congruent 
modulo AP to the Mth convergent of the continued fraction. Thus, it satisfies, 
modulo M^, a linear recurrence of order at most M. Hence it is eventually periodic 
modulo AI. 

(a) The estimate above of the order of the recurrence satisfied by modular reduc- 
tions of the pseudo- factorials can be vastly improved |8]. By the classical (2 x 2)- 
determinant identity of orthogonal polynomials and convergent denominators, the 
difference of two successive convergents of ( 32 1 satisfies the identity 

Pk+l{z) _ Pkjz) _ 0102 ■ ■ ■ ttkZ^'' 

Qk+i{z) Qk{z) Qk{z)Qk+i{z)' 

This specializes to the J-fraction (42 1 relative to F{z), when the aj are taken to 
be as in (43 1. By expressing that F{z) is the sum of the differences of its successive 
convergents, we then obtain, for any m > 0, 

Pm{.z) 0102 • • • flfcZ^*" 



{z) ^Qkiz)Qk+i{z)' 

In particular, since the aj are all integers and the Qj are integral with Qj{0) — 1, 
we have, 

P (z) 

F{z) = YTT^ *^)' = "i"2 ■ • • a™, 

nim \Z} 

in the sense that coefficients of both series are equal, after reduction modulo M . 
Thus, modulo M, the a„ satisfy a linear recurrence whose characteristic polynomial 
is exactly the denominator polynomial Qm{z) (reduced mod M). □ 

As an illustration, corresponding to to = 1, 2, 3, we find the congruences 

F{z) = -^ (mod3-l!2), F{z) = ,^—^ (mod 3 • 2!2), 
1 + z 1 + 2z^ 

n^)^ ,^'t_l'tfin3 (niod3^.3P), 
1 -|- 3z -I- K>z^ + lOz^ 

which already justify the data of Figure [3]for moduli 2, 3, 4, 6, 9, 12, 18; for instance, 
from the convergent Pi/Qi, we find = (—1)" (mod3). For m = 7, the form 
P7{z) ^ 3 + 6z + 5z2 -I- 2z3 + + 2z5 

rn^^^ ^"^^^'^ 

explains the observed patterns of («„) modulo 7 and the period equal to 36. By 
contrast, for to = 11, we find that 

P (z\ 

; = 1 + lOz + + 2z^ + + 4z5 + lOz*^ + 6z^ + + + z^" (mod 11) 
Qii(2:) 

(Qii reduces to 1 modulo 11), thereby establishing that the a„ with n > 11 are all 
divisible by 11. 

As the previous discussion suggests, congruence properties of pseudo-factorials 
are tightly linked to arithmetic properties of the Q polynomials whose exponen- 
tial generating function has been determined in Theorem [6] Let again 11^ denote 
an unspecified polynomial of degree r. Without attempting a general discussion, 
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we only remark here the existence of striking regularities, as summarized by the 
following data. First, for m a prime of the form 6/i + 1: 





mod 7 


mod 13 


mod 19 


mod 31 


P7 


He 




P19 _ Ills 


-P31 _ 1130 




1 + 4z6 


Ql3 " 1 + 11^12 


Qi9 " 1 + llzis 


Q31 ~ 1 + 4z30 • 



Finally, for m a prime of the form 6/i + 5: 



mod 


5 


mod 


11 


mod 


17 


mod 


23 


Q5 ~ 


n4 


Pii _ 

Qii ~ 




Pl7 _ 

Qn ~ 


ni6 


^23 _ 
Q23 " 


n22. 



10. Conclusion 



The relation between elliptic functions and continued fractions is an old subject, 
one that is especially rich. Connections are manifest with the theta function frame- 
work in the form of various types of q-series expansions, starting with Eisenstein 
and including the celebrated Rogers-Ramanujan identities [H Ch. 16]. Closer to 
our perspective are early contributions due to Stieltjes and Rogers regarding the 
Jacobian sn,cn framework. Conrad, first in his dissertation f5] then in collabora- 
tion with Flajolet [B], has elicited new connections with the Dixonian framework of 
the sm,cm functions. As the present work supplemented by further investigations 
of ours indicate, there are new continued fractions to be explored, attached to the 
WeierstraB and Dixonian frameworks. In this vein, we have recently discovered new 
elliptic continued fractions, relative to "equiharmonic" and "lemniscatic" numbers, 
which are lattice analogues of the pseudo-factorials — we plan to report on these in 
a future publication. 
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